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SECTION  1 
INTRODUCTION 


In  the  subject  area  of  analytical  process  modeling  for 
improved  composites  (APIC) ,  extensive  work  has  been  reported  in 
[1-4].  The  objective  of  the  APIC  program  was  to  develop  analytical 
capabilities  establishing  relationships  among  the  process  para¬ 
meters,  composite  properties  and  material  response.  This  capability 
can  provide  the  means  for  avoidance  of  fabrication  failures, 
processing  the  material  with  improved,  more  uniform  and  repeatable 
composite  properties,  and  help  in  further  understanding  and  modi¬ 
fying  the  process  sequence/cycle.  The  investigation  of  the  in¬ 
fluence  of  process  variables  on  composite  microstructure  may  lead 
to  rapid  and  low-cost  processes  which  produce  carbon-carbon 
materials  with  controlled  microstructure  and  desirable  properties. 

To  achieve  these  objectives,  two  models,  namely  the  process  en¬ 
vironment  model  (PEM)  and  the  mechanical  model  (MIPAC)  have  been 
developed  at  the  Science  Application  Inc.,  (SAI) ,  Irvine,  CA  [1]. 

The  PEM  deals  with  the  establishment  of  a  relationship  between 
the  physical  and  chemical  material  properties  of  constituent  materials 
and  processing  variables  (temperature,  temperature  rate  and 
pressure) ,  while  MIPAC  deals  with  the  prediction  of  the  effective 
thermoelastic  properties  and  strength  characteristics  of  the  final 
state  composite  material.  Computer  codes  have  been  developed  by 
SAI  for  both  PEM  and  MIPAC  models  for  the  in-process  description 
of  3-D  orthogonal  reinforced  composites.  For  this  report  we  have 
utilized  the  computer  code  version  available  in  October  1983. 

The  primary  purpose  of  the  PEM  is  to  predict  the  pressure 
and  temperature  environment  for  the  mechanical  model,  as  shown  in 
Figures  ia-b.  In  order  to  carry  out  such  an  analysis,  it  is 
necessary  to  account  for  the  process  thermochemistry,  heat  transfer, 
and  flow  field  for  the  liquid  and  gas  in  all  regions  of  the  process 
can.  Given  the  thermomechanical  environment,  the  mechanical  model 
then  oredicts  the  stresses,  strains,  and  displacements  in  the 
material,  followed  by  a  failure  analysis  to  determine  the  nature 
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II.  MESH  GENERATION 

A.  Mesh  Control  Card 


Format  (615) 


Columns 

Parameter  Description 

1-5 

MAXI 

maximum  value  of  I  in 
mesh*  MAXI  <  25 

6-10 

MAXJ 

maximum  value  of  J  in 
mesh*  MAXJ  £  100 

11-15 

number  of  line  segment  cards 

16-20 

ignored 

21-25 

number  of  material  block  cards 

26-30 

maximum  relaxation  iterations, 
if  blank  default  of  100  used 

*array 

dimensions  for  mesh  data 

Line 

Segment  Cards 

Format  (2(213,  2F8.3),  22X,  15) 

Columns 

Description 

1-3 

4-6 

7-14 

15-22 

23-25 

26-28 

29-36 

37-44 

67-71 


coordinate  of  first  point 

coordinate  of  first  point 
*  * 

(X)  coordinate  of  first  point 


coordinate  of  first  point 
coordinate  of  second  point 
coordinate  of  second  point 
( X)  +  coordinate  of  second  point 
coordinate  of  second  point 
line  segment  type  (0,  1  or  2) 

0  point  (first  point  only) 

1  connecting  straight  line  (vertical 
or  horizontal) 

2  connecting  straight  line  (+45°  to 
horizontal) 


in  Cartesian  coordinates 
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D. 


Time  Step  and  Output  Options 


Format  (4F10.0,  415) 


Columns 

Parameter 

Description 

1-10 

11-21 

TSTART 

TSTOP 

start  time  minutes  for 

_  ,  .  impregnation  & 

stop  time  u  .c  , 

r  hours  for  carb 

and  graph. 

21-30 

DT0 

time  step  interval 

31-40 

TPRINT 

minimum  print  interval 
selected  by  time 

41-45 

NPRINT 

minimum  print  interval 
selected  by  number  of  time 
steps 

46-50 

NMAXDT 

maximum  time  steps 

51-55 

NSK 

maximum  core  storage  (zero 
yields  default  of  9700) 

56-60 

NTTDOT 

0  no  nodal  output 

1  print  nodal  variables  only 

2  print  nodal  variables  and  time 

derivatives 

Equation 

Solution  Sequence  Options  Format  (415) 

Columns 

Description 

1-5 

number  of  equations  to  be  solved;  if 
blank,  sequence  1,  2  and  3 

6-10 

code 

for  first  equation 

11-15 

code 

for  second  equation 

16-20 

code 

for  third  equation 

codes 

;  1  temperature 

2  pressure 

3  gas  volume  fraction 

14 


1. 


ANALYSIS  CONTROL 


A. 

Number  of  cases  to  be  analyzed 

Format 

(15) 

B. 

Problem  Title 

Format 

( 8A10 ) 

C. 

Analysis  Options 

Format 

(915) 

Columns 

1-5 


6-10 


11-15 

16-20 


21-25 


26-30 


31-35 

36-40 

41-45 


Parameter 

NPP  0 

1 

I START  0 

1 
2 


Description 

axisymmetric 

plane 

start  at  beginning 
mesh  start  file* 
continuation  run** 


*MCC,  LS  and  MBA  cards  omitted 
**initial  value  card  also  omitted 


ISTOP 

IOPTN 

I  PROP 

I  QUAD 

NBF 

NTT 


0  run  complete  job 
1  generate  mesh  only 

1  carbonization 

2  impregnation 

3  graphitization 

0  material  properties  are 
not  temperature  dependent 

1  material  properties  are 
temperature  dependent 

0  material  properties  constant 
over  each  triangle 

1  material  properties  constant 
over  each  quadrilateral 

number  of  boundary  function 
cards 

number  of  time  tables 


NUMMAT 


number  of  materials 
(impregnant  and  preform) 
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APPENDIX  A 


PEM  COMPUTER  CODE  USER  INSTRUCTIONS 
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Reaction  midpoint 
temperature  and  its 
time  width. 
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where 

H  =  heat  of  reaction 

m^  =  rate  of  reaction/unit  initial  mass. 
Aw  =  total  weight  loss/unit  initial  mass. 

=  reaction  midpoint  temperature 
a  =  width  of  reaction  time 
s  =  gas  volume  fraction 


<D  = 

p2.  = 

f  (u) 

7. 


porosity 
liquid  density 

-  -±-  e'u2/2 

v'/2tT 

Thermophysical  properties 


of 


each  material: 


(a) 

(b) 

(c) 


Material  Porosity  (<£)•  ~ 

^  P 

Change  of  pitch  density  with  pressure  (^-)  • 

Material  property  table  for  each  temperature 
the  following  properties: 


(i) 

Thermal  conductivity  (K) 

(ii) 

Gas  molecular  weight  (Mq) 

(iii) 

Density  of  solid,  liquid  and  gas 

(vp 

(iv) 

Viscosity  of  liquid  and  gas  (Mj,, 

M  ) 
q 

(v) 

Specific  heats  of  solid, liquid  and  gas 

(vi) 

Permeability  (K^) 

• 

(vii) 

Material  emissivity. 

for 


£'pg} 


[C  C  0 
ps,  pf. 


,C 


pg 


) 
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.  .  O 


4.  Capillary  action  equation  constants 


a^  (i-l,...,4)  and  c  for 


P 


c 


c 


a  s  a 
1  1 


s  =  1 


-  s 


w 


s  =  wetting  fluid  saturation 
w  ^ 

(saturation  of  the  void  volume  with 
wetting  fluid) 


5.  Various  constants  required: 

(a)  Universal  gas  constant,  k  (.perfect  gas  law). 

(b)  Stefan-Boltzman  constant  (for  graph,  only). 

(c)  Gravitational  constant,  g. 

(d)  Energy  conversion  constant,  J  . 

(e)  Carbonization  rate,  T  . 

6.  Reaction  details  for  the  chemical  dissipation  function  (p^. 


* 

c 


3 

EH  .  m  .  ( 1  -  s )  $  p  , 

r  l  g  l  > 

i  =  l 


(for  carbonization  only) 


11,  pitch  decomposition 
i  =  | 2,  mesophase  formation 
[3,  carbon  formation 


SECTION  II 


MATERIAL  CONSTITUENT  PROPERTIES, 
ASSUMPTIONS  AND  OTHER  REQUIREMENTS  FOR  PEM 


A.  FIELD  THEORIES 

For  the  development  of  the  computer  code  the  following 
field  theories  are  used: 

a)  Fourier's  law  of  heat  conduction- 

fa)  Darcy's  law  of  flow  of  fluids  through  porous  media. 

For  further  details  the  reader  is  advised  to  consult  reference  [1]  . 

B.  INPUT  DATA  REQUIREMENTS  FOR  IMPREGNATION  ANALYSIS 

1.  Billet  (preform)  Dimensions 

2.  The  following  material  properties 

(a)  Change  of  pitch  density  with  pressure 


3p 

?P  ' 


(C 


i  Ipn 
p  3P 


C. 


(b)  Solid  material  porosity  (<t>) 


=  1 


^solid 

Vbulk 


(c)  Liquid  pitch  density  (p) 

(d)  Liquid  Pitch  viscosity  (y) 

(e)  Anisotropic  permeability  of  the  material  (KP) 

(f)  Boundary  conditions 

INPUT  DATA  REQUIREMENTS  FOR  CARBONIZATION  AND  GRAPHITI ZATION 
ANALYSES 


1.  Boundaries  of  billet  dimensions,  liquid  &  preform  constituents. 

2.  Process  cycle:  Distribution  of  temperature  &  pressure  with 
time . 

3.  Initial  and  boundary  conditions  for  temperature,  pressure 
and  gas  volume  fraction. 
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3 HR  HOLD 
AT  2300°C 


(Oo) 3HfUVa3dW31 


Figure  2b.  Graphitization  Schedule  for  APIC  Ideal 


TEMPERATURE 


Figure  2a.  APIC  Carbonization  Cycle. 


Processing  Sequence 


Weave  preform:  Carbon  fibers  are  stiffened  and  woven 
in  a  3-D  configuration. 

Cage  preform:  The  preform  is  placed  in  a  cage  to  allow 
adequate  access  of  the  impregnant  and  constrained 
gently  to  avoid  distortion. 

Preform  is  loaded  into  the  processing  can  and  fillers 
are  added. 

Impregnation:  Can  is  heated  to  a  prescribed  temperature 

in  the  impregnation  unit  and  the  impregnant  is  allowed 
to  flow  into  the  preform  at  temperature  and  pressure 
conditions . 

Carbonization:  The  can  with  impregnated  billet  is 

sealed/vented  and  subjected  to  a  prescribed  process 
cycle.  A  carbonization  process  cycle  used  by  SAI  [1] 
is  given  in  Figure  2a. 

Graphitization :  The  carbonized  billet  is  then  sub¬ 

jected  to  a  graphitization  process  cycle.  A  typical 
graphitization  process  cycle  is  given  in  Figure  2b. 

Last  three  steps  (step  4  through  step  6)  are  repeated 
until  desired  properties  are  achieved. 


and  extent  of  any  damage  which  might  occur.  This  must  be  repeated 
for  incrementally  applied  loading  until  the  entire  process  is 
modeled.  In  principle,  it  is  also  possible  that  the  output  of  the 
mechanical  model  can  alter  the  true  response  of  the  process  en¬ 
vironment  model  by  changing  material  thermal  conductivities  and 
permeabilities  at  local  damage  sites.  This  is  shown  by  the  feed¬ 
back  path  indicated  in  the  Figure  lb,  and  while  possible,  it  is 
not  a  real  time  link  between  computer  codes.  Approximate  satis¬ 
faction  of  this  link  can  be  provided,  however,  through  repetition 
of  analysis  and  the  use  of  engineering  judgement. 

In  this  report  we  will  concentrate  on  the  process  environ¬ 
ment  model.  All  the  required  (input)  material  properties  in  three 
different  phases  of  processing  i.e.  impregnation,  carbonization 
and  graphitization  are  listed  in  the  following  section.  User's 
instructions  for  the  PEM  computer  code  are  given  in  Appendix  A. 

These  instructions  are  for  the  analysis  of  impregnation, 
carbonization  and  graphitization  as  well  as  the  graphics  post  pro¬ 
cessor.  The  user's  instructions  for  computer  codes  for  3D  thermal 
and  diffusion  property  predictions  and  the  degree  of  graphitization 
analysis  are  not  given  in  this  report.  A  number  of  illustrative 
examples  are  treated  and  input/output  for  these  problems  are  ex¬ 
plained.  The  temperature  dependent  material  properties  used  in 
these  problems  are  exhibited  in  figures.  The  results  for  the  sample 
problems  are  presented  in  graphical  form.  Independently,  we  have 
derived  closed  form  solutions  for  impregnation  of  the  billet  in 
Cartesian  and  cylindrical  polar  coordinates  for  certain  special 
cases.  A  comparison  between  the  PEM  results  and  closed  form  solution 
results  has  been  made.  There  exists  a  good  agreement  between  the 
converged  results  obtained  by  the  PEM  and  closed  form  solution  for 
both  the  rectangular  and  cylindrical  polar  coordinate  systems. 

For  the  complete  processing  of  3-D  carbon-carbon  material 
the  following  processing  sequence  may  be  employed: 


CONTROLLABLE  PROCESS 


.  •  V  • 


c  LxJ 

,  —  2 

02  o 

Cd 

rv-  — - 


^  o 

LU 

2  O  uj 
O  —  H- 

_  <C  i — 

>  s 1/1 

g  o;  lu 
LU  D_  O 
O  *< 

2:  lu  < 

E  |  < 


oo  o 

LU  2! 

2  => 

p*  q;  2 

H  b  ^ 

^  <r 


< 
etc  o 


<d.  ~T~  <d  C_J 
<  O  O 

tt:  o  —  2: 
i-  5  on 
.  <  DQ 

h-  a;  u_  — 

<:  lu  u_  ~3 

n:  — 

m  t —  q  lj_ 

•  •  •  • 


LU  ^ 

oo  O 

^  ^  LO 
O  CQ  | — 
H  Q  D 
l/l  LU  CL 
i  i  i  LU  — ^ 

a!  hb  _ 


LU 

Cd 

LU 

Qt 

ZD 

1 — 
< 

ZD 

Cd 

LO 

LU 

OO 

LO 

Q_ 

UJ 

LU 

2 

1— 

Cd 

UJ 

<C  i 

Cl. 

1— 

Cd 

• 

• 

• 

t-  2 

y  e= 

— J  o 
—  c 
00  at 
b— 
2  oo 
O  CQ 

on  <' 
~  u_ 

5  ° 


3 


Figure  lb.  Relationship  Between  Process  Environment  Model  and 
Mechanical  Model. 


FINITE  ELEMENT  TRANSIENT  FIELD  SOLUTIONS 
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DEGREE  OF  GRAPH  IT  I ZATI  ON  ANALYSIS 


C.  Material  Block  Assignment  Cards  Format  (515,  F10.0) 

(material  outside  billet  must  be  material  number  one) 

Columns  Description 


1-5 

material  definition 

number  (1-12) 

6-10 

minimum  I  boundary 

11-15 

maximum  I  boundary 

16-20 

minimum  J  boundary 

21-y5 

maximum  J  boundary 

26-35 

material  property  inclination  angle 
in  degrees  (usually  blank) 

Note:  the  orientation  of  the  M-N-P 

material  coordinate  system  with 
respect  to  the  R-Z-P  body  coordinate 
system  is  shown  below.  The  angle 
is  input  via  a  MATERIAL  BLOCK 
ASSIGNMENT  CARD.  Both  are  right  handed 
systems . 


Z 


BOUNDARY  FUNCTIONS 


Format  (615,  2F10.0) 


Each  card  defines  input  data  for  a  group  of  nodal  points 
bounded  by  1 1 ,  12,  Jl,  J2,  on  the  external  boundary  of  the 
body.  For  a  line,  II  =  12  or  Jl  =  J2 .  For  a  point,  II  = 


12  and  Jl  =  J2.  The  time  variation  for  each  card  input 

below  is  specified 

by  TIME  TABLE  CARDS. 

Columns 

Description 

1-5 

KIND,  the  type  of  boundary  function 
(1  -  temp.,  2  -  pressure,  3  -  gas 
volume  fraction  4  -  graphitization 
source  temp.) 

6-10 

11,  the  initial  value  of  I 

11-15 

Jl,  the  initial  value  of  J 

16-20 

12,  the  final  value  of  I 

21-25 

J2,  the  final  value  of  J 

26-30 

NTABLV ,  corresponding  time  table 
identification  number  (zero  if 
function  does  not  vary  with  time) 

31-40 

VI,  the  value  of  the  function  at 
(11,  Jl) 

41-50 

V2,  the  value  of  the  function  at 
(12,  J2 ) 

INITIAL  VALUES 

Format  (30X,  3F10.0) 

The  initial  value  of  temperature,  pressure,  and/or 
gas  volume  fraction  for  each  equation  to  be  solved, 
input  in  problem  sequence  order. 

Columns  Description 


31-40 

41-50 

51-60 


constant  initial  value,  first  equation 
constant  initial  value,  second  equation 
constant  initial  value,  third  equation 


AUXILIARY  EQUATION  CONSTANTS 


A.  Capillary  Action  Equation  Format  (5F10.0) 


Columns 

Descri 

ption 

1-10 

capillary 

pressure 

constant , 

al 

11-20 

capillary 

pressure 

constant , 

a2 

21-30 

capillary 

pressure 

constant , 

a3 

31-40 

capillary 

pressure 

constant. 

a4 

41-50 

capillary 

pressure 

strength , 

c 

Pc  »  c  ai  s  +  a2t(i~=^)  ‘  a4] 


(a-^  to  a^  are  dimensionless  and  c  is  in  PSI) 


B.  Material  Constants 


Format  (7F10.0) 


Columns 


1-10 

11-20 

21-30 

31-40 


41-50 

51-60 

61-70 

71-80 


Description 

universal  gas  constant,  R  -(BTU/lb.  mole-°R) 

residual  liquid,  S*r  (normally  blank) 

residual  gas,  S  (normally  blank) 

9r  2  4 

Stefan-Boltzman  constant  -  (BTU/Hr  in  °R  ) 

(graphitization  problems  only) 
gravitional  acceleration,  g  -  (in/Hr^) 
gravitional  constant,  gc  -  (in/lb-F  Hr  ) 
energy  conversion  constant,  Jc(in  lb-f/BTU) 
c a rboni nation  rate.  T 


MATERIAL  THERMOPHYSICAL  PROPERTIES 
(one  set  of  input  properties  per  material) 


A.  Material 

Columns 

1-5 

6-10 

11-15 


16-20 


21-30 

31-40 

41-50 


Identification  Format  (415,  3F10.0) 

Description 

identification  number  (1-12) 
number  of  temperatures  (1-12) 
isotropy  parameter 

0  anisotropic  material 
1  isotropic  material 
material  type 

1  solid 

2  liquid 

3  porous  (normally  used) 
change  of  pitch  density  with  pressure 
material  porosity 

solid  material  density 


Material  Property  Table 


Format  (7F10.0/10X,  5F10.0/ 
10X,  6F10.0) 


11-20 

21-30 

31-40 

41-50 

51-60 


gas  viscosity, 

liquid  pitch  viscosity,  y^ 

yarn  or  filler  specific  heat,  Cps 

permeability,  or 

permeability, 


Th i rd  Card 


11-20 

gas  density, 

21-30 

gas  specific  heat,  C 

pg 

31-40 

liquid  pitch  specific  heat,  C  „ 

P  ° 

41-50 

material  emissivity, 

Units  for 
8-19  also 

these  quantities  are  mentioned  in  Figures 
see  Tables  3  and  4. 

VII. 


TIME  TABLE  INFORMATION 


A.  Time  Table  Identification 

Columns 

1-5 
6-10 


B.  Time  Tables 
Columns 

1-10 

11-20 

21-30 

31-40 

VIII.  TERMINATION  OF  INPUT 

A.  End  of  Case  Card 
enter,  "END  OF  CASE 

B.  End  of  Data  Card 
enter, "END  OF  DATA" 


Format  (215) 

Description 

identification  number 
number  of  ordered  pairs  (<_  20) 

Format  (8F10.0) 

Description 

first  table  time 
first  function  value 
second  table  time 
second  function  value 

Format  (2A10) 


Format  (2A10) 


APPENDIX  B 


POST-PROCESSOR  AND  GRAPHICS  CODE  USER  INSTRUCTIONS 


Post  processor  and  graphics  code  developed  by  SAI  [1] 
is  helpful  in  obtaining  contour  plots  of  various  parameters 
during  processing.  This  code  is  written  such  that  a  SC-4020 
plotter  is  used  -  and  therefore  requires  the  relevant  plot  routine 
library  attached.  Data  preparation  instructions  are  as  follows: 

1.  TITLE  AND  CONTROL 


A.  Title 

Format 

( 8A10 ) 

B.  Options  Control 

Format 

(515) 

Columns 


Description 


1-5 


IDEN 


0  no  density  calculation 
1  calculate  density 


6-10 


IMESH 


0  no  mesh  plot 
1  plot  mesh 


11-15 


ICON 

0 

no  contour  plot 

1 

plot  contours 

16-20 


I  PLOT 


0  standard  plots 
1  non-standard  plots 

(plot  scale  card  required) 


21-25 


I START  0 
1 


data  on  tape  8 

data  on  tapes  21,22,23,  and  24 
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II. 


MATERIAL  PROPERTY  Format  (2I5/2F10.0) 

(omitted  if  density  is  not  calculated) 


A.  Material  Identification 


Columns  Description 

First  Card 

1-5  number  of  materials 

6-10  number  of  ordered  pairs  in  pitch 

density  table  (<_  12) 

Remaining  Cards  -  (one  card  each  material) 
1-10  porosity 

11-20  solid  constituent  density 

B.  Pitch  Density  Table  Format  (8F10.0) 


Columns 


Description 


1-10 

11-20 

21-30 

31-40 


first  table  temperature 
first  table  density 
second  table  temperature 
second  table  density 
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PLOT  SCALE  OPTION  Format  (5F10.0) 

(omitted  if  standard  plot  selected) 


Columns 

Description 

1-10 

minimum 

R 

coordinate 

plotted 

11-20 

maximum 

R 

coordinate 

plotted 

21-30 

minimum 

Z 

coordinate 

plotted 

31-40 

maximum 

Z 

coordinate 

plotted 

41-50 

rotation  ; 

parameter 

(if  1.0,  Z-axis  is  vertical; 
if  0.0,  Z-axis  is  horizontal) 

CONTOUR  PLOTTING  OPTION 

(omitted  if  contour  plots  not  selected) 

A.  Plot  Control  Format  (315) 


Plot  Function  Options 

One  card  is  required  for  each  variable  to  be  plotted. 
A  card  may  utilize  either  of  the  two  options  defined 
below. 

i.  Automatic  Contour  Selection  Format  (815,  2A10) 

Columns  Description 


1-5 

plot  code  1  temperature 

2  pressure 

3  liquid  fraction 

6-10 

number  of  contours  desired 
(enter  as  negative.  <  10) 

11-15 

first  of  three  acceptable  contour 
divisions  (default =  1) 

16-20 

second  of  three  acceptable  contour 
divisions  (default =  2) 

21-25 

last  of  three  acceptable  contour 
divisions  (default  =  5) 

26-30 

variable  range  plotted 

0  entire  range 
1  non-negative  range 
-1  non-positive  range 

31-35  symbol  frequency  (default  =  4) 

36-40  range  selection  method 

0  based  on  current  time 

1  based  on  total  time 

(Tape  8  must  be  used,  not 
Tapes  21  ,  22,  23  and  24  I START 

variable's  units,  printed  on  plot 

time  units,  printed  on  plot  (Hours) 


41-50 

51-60 


Specified  Contour  Values  Format  ( 215 , 30X ,  2A1 0/8F10 . 0 ) 


Columns 


Description 


6-10 

41-50 

51-60 


First  Card 


plot  code 


temperature 

pressure 


3  liquid  volume  fraction 

4  density 

number  of  contours  (<_  10) 
variable's  units,  printed  on  plot 
time  units,  printed  on  plot 


Second  Card 


1-10 

11-20 


first  contour  to  be  plotted 
second  contour  to  be  plotted 


TERMINATION  OF  INPUT 


End  of  Data  Card 


Format  (2A10)  enter  "END  OF  DATA" 
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APPENDIX  C 
PEM  ILLUSTRATIONS 
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F  i  'h:  t'o  9.  Liquid  Thermal  Conductivity  KNN 
Temperature  T  (°F). 


Versus 


I 


T(°F) 


i 


Figure  8.  Liquid  Thermal  Conductivity  K,„, 
'  1  MM 

Temperature  T  (°F). 
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fence  F 

- 

0  2 

. b B COE -02 

.  1500E+030 

4fe  505 

- 

01- 

1  0  045 

- 

G4 

4  0 

Efei 

- 

1  ? 

.<-0  135-04 

.29805+00 

10  0Gr 

- 

0 

r 

D 

13  0  3  5 

- 

06 

fe  L 

3,  H  l 

- 

C4 

.20334+30 

.32505+000 

h  6 

:  J  r 

- 

fe  <ft 

.  F  ft  0  C  _  -  G  2 

.15005+330 

4d  3  05 

- 

0  1- 

1  0  C  4  fe 

- 

04 

71 

4  7l 

- 

4_ 

.:C13E-04 

. 31235+03 

moor 

- 

0 

D 

1C  C35 

- 

06 

L.  ‘ 

7  fe  L 

- 

•  fe  1  t  OL  ♦  0  0 

•38535+000 

fe  •• 

r  o  - 

- 

'  ' 

.  -  fe.  0  3  fe  -  0  2 

.15005+030 

44  73  5 

- 

01- 

55305 

- 

05 

* 

'  >  ?  •" 

- 

13 

.fe'-tie.fe-CE 

.52235+00 

1  0  3  0  4 

- 

n 

O 

1  0  0  3  5 

- 

0  6 

-  t 

7  7  fe 

- 

fe  4 

.  2  2  4  C  4  ♦  C  0 

.12175+010 

fe  - 

fe  0  t 

- 

Z  xi 

.  fe  8  C  3  i  -  G  2 

.  15  005  +  G33 

44  4  GE 

- 

C  1- 

719  0  5 

- 

Ofc 

hr 

1  /  L 

- 

i : 

•  25365-0  8 

•  33G0t+uLi 

10  01  E 

- 

0  6 

10  LD5 

- 

0b 

4  1 

'  fe  5 

- 

J4 

. 2 35  OF  +  C 0 

,4,<hO5  +  C00 

30 

COE 

- 

G  1 

.  30  0  05-01 

.  150  C5  +  C  30 

44  105 

- 

01 

1  0  7  0  5 

- 

04 

1  U 

35  4 

- 

11 

.44<l9E-G7 

.  34b  e.  5  +  00 

5^  02r 

- 

0 

P 

50  03  5 

- 

08 

33 

794 

- 

04 

.  1 4  8  9  5  +  0  0 

.  39265  +  C  3  0 

4  0  0  0  5 

- 

0  1 

.  4(?  0CE  -  0  1 

•  13  6  05+020 

4  3  9  0  5 

- 

Cl 

2  3  c  C  5 

- 

04 

1  1 

n  f,  •: 

- 

1  1 

. _C 1 35- 04 

•36604+00 

1  fe  C  0  fe 

- 

n 

-j 

1  0  0  0  5 

- 

C  9 

2  u 

F-  0  E 

- 

04 

. ,  6  2.  3  E  +  C  0 

.51705+000 

5  b 

6  0  fe 

- 

C  1 

.  5  -fe  G  3  -  C  1 

.75005+010 

4  3  655 

- 

0  1 

253  0  5 

- 

04 

12 

fed 

- 

1 1 

.  «.  C:  1  3  4  -  0  4 

.37505+30 

10  COE 

- 

3  9 

10  COE 

- 

09 

t  fe 

9  1  fe 

- 

34 

.  4  6  5  3  i  ♦  0  0 

.50505+000 

5  6 

COL 

- 

fe  1 

.  -  4fe  0L- 31 

.22805+020 

4  3  5  3  5 

- 

01 

fe  6  3  0  5 

- 

06 

1  5 

.  Q  r 

- 

1 1 

.20135-04 

•  384  OF +  00 

1  0  0  3  5 

- 

r 

9 

1  D  0  3  5 

- 

09 

2  L; 

2  1  fe 

- 

0  4 

•  2  6  i  2  t  ♦  G  0 

.47605+300 

.  •  b  C  'J  t 

- 

.  1 

.  5  fc  0  0 1  -  0  1 

•10705+050 

4  3  3  5  5 

- 

0  1- 

655  OF 

- 

05 

1  3 

■  *8l 

- 

fe  1 

.  _  3 1 3L -  34 

.39204+00 

1  0  00  5 

- 

0 

i 

1  3  0  0  5 

- 

0  9 

fe  l 

fe  fe 

- 

fe  4 

. 971Cfe+GC 

•430 05+000 

i  <L' 

fe  0  L 

- 

'J  1 

.72335-31 

.47705  +  1/20 

4  3  20  5 

- 

01- 

86  DC  fe 

- 

0  5 

i  H 

4  9  fe 

- 

1  1 

.  0  1  3  L  -  G  4 

.40005+00 

1  C  0  G  L 

- 

n 

o 

- 

1  0  0  3  5 

- 

3  9 

n  . 

i  :l 

- 

0  4 

.  2  7>T-  *00 

.4jCuE  +  0  00 

a  r 

.'J  C  fe 

- 

0  1 

. -0CG5-G1 

•27805+320 

4  3  3  j  E 

- 

01  - 

9  6*  0  F* 

- 

05 

i  Z> 

1  4  fe 

- 

i  : 

.<-01  '5-04 

.40775+00 

roL 

- 

** 

1  0  C  C  5 

- 

r  ^ 

sz 

9  OF 

- 

0  4 

. ;;  7  7  c  ♦  3  c 

.4^/774  +  uo0 

i  4 

0  c 

♦ 

3  r 

.  1  4  ?  fe  5  ♦  0  0 

.16304+020 

6  4  3  3  5 

- 

0  1 

2  ,  -  0  c 

- 

04 

i  r, 

fe  4  F 

- 

i  1 

. ..  c:  3-:-G4 

.45335+00 

1  0  0E  L 

- 

0  9 

1  0  00  5 

- 

09 

t 

2  Mfe 

- 

0  4 

.  :  g  4 1  -:  ♦  o  o 

.45305+300 

4  2 


IX,  Z)  (0,0)  (2,0) 

II,  J)  (1,1)  (5,1) 


Figure  7. 
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Model . 
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and  permeability  (Table  5) .  The  contour  plots  of  all  the 
significant  parameters  at  different  times  during  processing  for 
the  three  cases  are  given  in  Figures  26-58.  For  each  case, 
temperature  contours  are  given  at  six  values  of  time,  starting 
from  T  =  5  hours  to  T  =  30  hours, with  a  gap  of  five  hours.  Pressure 
and  density  contours  are  plotted  at  two  times,  T  =  5  hours  and 
T  -  30  hours.  Also,  the  temporal  variation  in  temperature, 
pressure,  density  and  gas  volume  fraction  at  point  (4", 6")  has 
been  shown  for  each  case  (Figures  36,  47,  58).  Those  are  not 
a  part  of  the  TEM  code.  T’-e  carbonization  and  rraphiv.izat ion 
problems  solved  here  are  independent  of  each  other. 
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P  ( 2  ,  z ,  t )  =  90  PSI  0  <  z<  10" 
P  (x, 1 0 V t )  =  90  PSI  0  <  x  <  2" 


i£  |  =  o 

3z  Z  =  0 
p(x ,  z  ,  0  )  =  0 


The  input  data  for  this  analysis  is  given  in  Table  2.  The  pressure 
distributions  as  a  function  of  spacial  coordinates  and  time  are 
given  in  Figures  4  and  5. 

For  carbonization  analysis  two  impregnated  billets,  as 
shown  in  Figure  6,  XZ-plane,  are  considered  to  be  placed  in  a 
processing  can.  Because  of  geometric  symmetry,  half  of  the  can  was 
modeled.  The  boundary  conditions  imposed  in  this  problem  are  also 
shown  in  Figure  6.  The  finite  element  model  shown  in  Figure  7 
has  been  used.  The  temperature  dependent  material  properties  of 
the  billet  preform  and  liquid  pitch  are  given  in  Tables  3,  4  and 
5.  The  properties  given  in  Tables  3  and  4  are  given  in  graphical 
form  in  Figures  8-23.  The  carbonization  process  temperature 
schedule  is  given  in  Figure  24.  The  carbonization  analysis  is 
done  up  to  30  hours  of  processing.  Table  1.  Because  of  the  length 
of  the  processing  time  and  the  size  of  the  problem,  the  analysis 
was  conducted  in  two  runs  for  each  set  of  preform  porosity  and 
density  values,  i.e.  (1)  0  £  T  £  15  (hrs),  (2)  15  £  T  £  30  (hrs) 
as  discussed  earlier.  The  input  data  for  the  carbonization  pro¬ 
cessing  analysis  up  to  15  hours  is  given  in  Table  6.  The  data 
generated  in  that  run  was  stored  in  two  tapes,  SAVE  =  SAVE  1  and 
RESTRT  =  RSTRT  1.  These  two  tapes  are  called  to  be  used  in  the 
next  run  for  conducting  the  subsequent  processing  analysis  for 
15  ■_  T  _  3  0  (hours).  The  input  data  used  for  this  analysis  is  niven 
in  Table  7.  Similar  data  sets  were  made  for  carbonization  analysis 
using  other  two  sets  of  porosity  (.55,  .7),  density  (.065,  .04) 


TABLE  1 

PROBLEMS  SOLVED  FOR  ILLUSTRATING  THE  USE  OF  PEM  COMPUTER  CODE 


In  Appendix  A,  the  user  instructions  are  given  for  the  first 
four  models.  The  input  data  and  the  results  for  a  number  of 
examples  follow  . 

The  problems  solved  for  illustrating  the  use  of  the  PEM 
computer  code  are  given  in  Table  1.  The  preform  and  the  pro¬ 
cessing  can  dimensions  for  different  processes  are  given  in  Figures 
3,  6  and  7.  The  PEM  has  been  used  to  study  the  impregnation, 
the  carbonization  and  the  graphitization  analyses.  For  each 
case,  the  finite  element  model  and  the  relevant  boundary/initial 
conditions  are  described  separately.  The  thermochemical  and 
physical  properties  of  the  billet  preform  and  the  liquid  pitch 
are  provided  by  SAI  [5].  In  modeling  the  billet  two  dimensional 
Cartesian  coordinates,  XZ,  have  been  considered.  The  impregnation 
and  the  graphitization  analyses  are  conducted  for  one  set  of 
porous  medium  porosity  and  density  values.  The  carbonization 
analysis  is  done  for  three  cases  of  porosity  and  density  of  the 
preform.  The  computer  program  was  run  for  15  time  steps  each 
time  and  the  restart  capability  was  used  for  the  30  hour  carboni¬ 
zation  analysis.  In  the  first  run,  i.e.  0  £  time  <_  15  hours, 
the  output  was  stored  for  its  usage  in  the  second  run  which  is 
from  15  hours  to  30  hours.  All  the  results  obtained  during 
execution  of  these  programs  were  stored  for  obtaining  contour 
plots . 

In  the  finite  element  modeling  of  impregnation  process, 
Darcy's  law  for  the  flow  of  compressible  fluid  through  porous  media 
has  been  used.  The  effect  of  an  external  applied  pressure  has 
been  considered  and  that  of  the  capillary  pressure  has  been  ig¬ 
nored.  Also,  the  time  taken  in  pouring  the  fluid  in  the  impreg¬ 
nation  vessel  has  been  assumed  negligible.  Figure  3  shows  the 
billet  finite  element  grid  with  relevant  boundary  conditions. 
Because  of  the  geometric  symmetry  a  quarter  of  the  billet  has 
been  considered  for  finite  element  analysis.  Thus,  the  following 
boundary  and  initial  conditions  were  employed: 
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As  shown  in  Figure  la,  the  PEM  Code  has  the  following 
capabilities  [1]: 

1.  Impregnation:  This  code  is  based  upon  Darcy's  law  of 
flow  of  compressible  fluids  through  a  porous  media. 

The  input  parameters  required  for  this  analysis  are 
given  on  page  8.  Pressure  boundary  conditions  are 
assigned.  The  effect  of  mechanical  or  capillary  pressure- 
can  be  studied  separately.  The  present  model  does  not 
account  for  the  combined  effect  of  capillary  and 
mechanical  pressure  conditions. 


2.  Carbonization:  This  is  the  most  involved  process  in 
the  PEM  and  requires  a  unified  treatment  of  the  trans¬ 
fer  of  thermal  energy  by  considering  heat  transfer, 
fluid  motion  and  thermochemical  reactions.  Darcy's 
formulation  has  been  utilized  in  developing  this  part 
of  the  code. 


f 


3.  Graphitization :  This  capability  of  the  PEM  code  is 
based  upon  the  heat  conduction  in  the  billet  including 
the  influence  of  internal  heat  generation  due  to  the 
carbonization  and  graphitization  during  processing. 

4.  Postprocessor:  This  code  enables  the  user  to  obtain 
the  results  of  different  parameters  in  the  form  of 
contour  plots. 

5.  TEMP3:  This  code  has  the  capability  to  synthesize 
the  thermal  conductivity  matrix  [k]  and  the  diffusion 
matrix  [ D ]  for  a  3D  composite  given  the  similar  pro¬ 
perties  for  the  constituents  of  the  composite. 

6.  GRAPH:  The  degree  of  graphitization  of  the  material 
in  processing  at  any  point  of  time  can  be  estimated 
by  using  this  model. 
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Figure  11.  Liquid  Pitch  Density  p  f  Versus 
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Figure  17.  Gas  Density  p  Versus  Temperature  T. 
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Figure  20.  Solid  Thermal  Conductivity  K 
Temperature  T. 
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CARBONIZATION  ANALYSIS  OUTPUT 
PREFORM  POROSITY  =0.31 
PREFORM  DENSITY  =  0.1 
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In  the  input  of  Table  6,  the  processing  analysis  was  done 
up  to  15  hours.  The  data  generated  in  that  run  was  stored  in 
two  tapes,  SAVE  =  SAVE  1  and  RESTRT  =  RESTRT  1.  These  two  tapes 
are  attached  to  use  in  the  next  run  to  conduct  the  processing 
analysis  from  15  hours  to  30  hours.  The  input  for  this  analysis 
is  given  in  Table  7. 

The  aforementioned  data  sets  are  given  for  a  billet  with 
porosity  =  0.31,  density  =  0.1  and  corresponding  values  of 
permeability  from  Table  5.  Similar  data  sets  with  other  two 
values  of  porosity  (0.55,  0.7),  density  (0.065,  0.04)  and  per¬ 
meability  (Table  5)  were  used  to  conduct  carbonization  analysis. 
The  contour  plots  of  significant  parameters  for  all  the  three 
cases  are  given  in  the  following  section.  In  the  plotting  run 
for  these  contour  plots,  we  have  used  the  automatic  scale  choice 
option,  and  for  that  reason,  the  figures  are  not  consistent  with 
the  finite  element  grid  axis  shown  in  Figure  7.  Contour  figures 
are  consistent  with  the  billet  axis  and  boundary  conditions  shown 
in  Figure  25. 
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IMAGES  OF  DATA  FOR  CARBONIZATION  ANALYSIS  FOR  TIME  T=15  TO  T=30  HOURS 
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Cards  63  through  98  are  the  same  as  solid  preform  properties  Table  4  on  page  55. 


The  following  material  properties  are  common  to  both  (Table  3) 
liquid  pitch  material  property  and  (Table  4)  solid  preform  material 


property. 


Gas  molecular  weight  M  (Fig.  #10) 

g 

Liquid  pitch  density  (Fig.  #11) 

Rate  of  change  of  pitch  density  with  respect  to 


temperature  — ^7 
Gas  viscosity  y 


(Fig.  #12) 
(Fig.  #13) 


Liquid  pitch  viscosity  y^  (Fig.  14) 

Gas  density  p  .  (Fig.  #17) 

g 

Gas  specific  heat  C  (Fig.  #18) 

pg 

Liquid  pitch  specific  heat  C  .  (Fig.  #19) 


I  ’C) 


Fiquro  23.  Solid  Permeability  (or  K? 

Temperature.  x  z 
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PLOT  SAV5 
TEMPERATURE 


Figure  30.  Temperature  Contours  in  Carbonization  Analysis 


Figure  31.  Temperature  Contours  in  Carbonization  Analysis 
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Material  Density  Contours  in  Carbonization  An 


100 


250 


CARBONIZATION  ANALYSIS  OUTPUT 
PREFORM  POROSITY  =  0.55 
PREFORM  DENSITY  =  0.065 
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Analysis 


/////MV 


Carbonization  Analysis 
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.  J.. 


Temperature  Contours  in  Carbonization  Analysis. 


[’ornpuraturc  Contours  in  Carbonization  Analysis. 
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PRESSURE 


Material  Density  Contours  in  Carbonization  Analysis 
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CARBONIZATION  ANALYSIS  OUTPUT 
PREFORM  POROSITY  =0.7 
PREFORM  DENSITY  =  0.04 
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0.65 


T  ( °F) 


Figure  63. 


Solid/Filler  Specific  Heat  C 
Versus  Temperature  T.  ps 
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Figure  62.  Liquid  Pitch  Density  p„  Versus 
Temperature  T. 


Figure  61.  Solid  Thermal  Conductivity  Versus 

Temperature  T.  1 


T  ( °F) 


Figure  60.  Solid  Thermal  Conductivity  Versus 

Temperature  T. 


TABLE  9 


CARBONIZED  SOLID  MATERIAL  PROPERTY  TABLE  CORRESPONDING 
TO  THE  INPUT  DATA  QUANTITIES  GIVEN  ON  PAGE  21,  i.e. 


T’  lW  'W 

Bp. 

g*  3T 

Format 

(4F10.0,  10X, 

F10.0) 

53  7. 

c 

ps 

c  ,  c 

pg  p£ 

at  each  temperature  T  is: 

3.9688  3.6242 

Format 

Format 

16. 

(30X,  F10.0) 

(30X,  2F10.0) 

0 .0  5058 

1302. 

2.6231 

2.4380 

0.360 

0.500 

16. 

0.82 

0.05058 

1  £.6  2. 

2.3866 

2.228 

0.380 

0.500 

16. 

0 .82 

0 .05058 

2  292. 

2.2261 

2 .2256 

0.4  02 
0.500 
16. 

0.82 

3  .  0  6322 

2  460  . 

2.0937 

2 .0936 

0.451 

0.490 

16. 

0.82 

0.06503 

2660  . 

1.9284 

1 .9283 

0.462 

0.485 

16. 

Q.R2 

0  .0  7226 

3060. 

1. 6793 

1 .6793 

0.468 
0.4  68 
16. 

0  •  q  2 

0.07226 

346  0  . 

1.734 

1 . 734 

0.475 

0.475 

16. 

0.82 

0.07226 

4  46  0  . 

1 .5  76 

1.578 

C  .505 
0.505 

16. 

0.82 

0.07226 

546  0. 

1.494 

1.494 

0.521 

0.521 

16. 

•J  .82 

0.07226 

6  4  b  0  • 

1.464 

1 .464 

0.525 

0.525 

16. 

0.8  2 

0.07226 

0.525 

0.525 

0  ,°2 

In  the  following  Figures  the  temperature-dependent  parameters  from 
the  foregoing  table  are  given.  Graphitization  temperature  variation  with  time 
is  also  given.  Similar  plots  could  be  obtained  for  carbonized  liquid  material, 
since  in  the  case  considered  here,  the  difference  between  the  material  properties 
for  carbonized  liquid  and  carbonized  solid  is  small,  material  property  plots 

are  civen  only  for  carbonized  solid. 
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1  vt:- 


■}  . 

7  . 

19. 


TABLE  8 

CARBONIZED  LIQUID  MATERIAL  PROPERTY  TABLE  CORRESPONDING 
TO  THE  INPUT  DATA  QUANTITIES  GIVEN  ON  PAGE  21,  i.e. 


T, 

W  Sin’  ms’ 

8T 

Format 

(4F10.0,  10X, 

F10.0) 

C 

ps 

C 

Pg 

For 

537. 

’  cp£ 

each  value  of 

3.795 

T  are: 

3.711 

Format 

Format 

16. 

(30X,  F10.0) 

(30X,  2F10.0) 

0 .05058 

1  302  . 

2.507 

2.458 

0.360 

0.500 

16. 

0.8  2 

0.05058 

1  562. 

2.281 

2.238 

0.360 
0.50  0 
16. 

0.8  2 

0.05058 

2  292  . 

2  •  2  j.2 

2.211 

0.4Q2 

0.500 

16. 

0.  «2 

0.06322 

2460. 

2.0852 

2.0850 

0.451 

0.490 

1  6  • 

0.82 

0.06503 

2660. 

1.9245 

1.9245 

0.462 

0.485 

16. 

0.R2 

0.07226 

3060. 

1.6791 

1.6791 

0.468 

0.466 

16. 

0.R2 

0 . 0  7226 

3460  . 

1 .734 

1  .734 

0.475 

0.475 

16. 

0.82 

0.0  7226 

4460. 

1.578 

1  .576 

0.505 

C  .505 

16. 

0  •  p  2 

0 .0  7226 

5460. 

1.494 

1.494 

0.521 

0.521 

16. 

0.82 

0 .07226 

6460. 

1.464 

1.464 

0.525 

0.525 

16. 

0.«2 

0.07226 

Curing  time  table 

0.525 

0.525  0.82 

for  graphitization  is: 

Time,  Temp 

537.  1. 

Format 

1031.4 

(8F10.0) 

3. 

1  67  ?.  4 

5. 

2291.4 

11. 

2615.4 

13. 

2813.4 

16. 

7101.4 

2  3. 

3353.4 

29. 

5225.4 

3  0. 

In 

these  tables  T 

is  given 

in  °R  where 

as  in  plots  T 

is  given 

in 

(°F). 

2111.  4 

3011.4 

5387.4 
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For  graphitization  analysis,  the  same  two-billet  can  used 
for  carbonization.  Figure  6,  has  been  considered.  The  finite 
element  grid  and  the  boundary  conditions  shown  in  Figure  6  have 
been  used.  The  temperature  dependent  material  properties  of  the 
billet  preform  and  liquid  pitch  after  carbonization  are  given  in 
Tables  8  and  9.  These  properties  do  not  necessarily  correspond  to 
any  of  the  three  cases  of  billets  analyzed  for  carbonization  process 
in  this  report.  The  graphitization  process  temperature  schedule 
is  given  in  Table  8  and  also  in  Figure  65.  The  PEM  computer  input 
data  for  this  problem  is  given  in  Table  10.  The  processing 
analysis  is  conducted  up  to  30  hours.  Figure  66  shows  the 
temperature  profile  results  at  two  locations  of  the  billet  during 
graphitization . 
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Figure  58.  Variation  of  Response  Parameters  with  Time  at  the 
Center  of  the  Modeled  Billet  i.e.  Coordinates 
(4",  6")  . 
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Figure  56.  Material  Density  Contours  in  Carbonization  Analysis 


Figure  64. 


Liquid  Pitch  Specific  Heat  C  ^ 
Versus  Temperature  T.  ^ 
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Figure  65 


TIME  (HRS) 


Temperature  at  Different  Times  at 
Two  Locations  in  the  Billet  Durinq 
Graphitization. 


Figure  66. 


Governing  Equations,  Solutions  and  Results 


Darcy's  law  for  the  flow  of  a  compressible  fluid  through  a 
porous  media  has  been  used  for  the  impregnation  analysis.  The 
relevant  relations  [6]  for  this  investigation  are  presented  in 
this  section.  The  governing  equation  of  flow  of  a  compressible 
fluid  through  porous  media,  in  terms  of  pressure  in  the  medium, 
will  be  solved  for  the  case  of  uniform,  isotropic  and  time- 
independent  properties.  The  following  three  types  of  boundary 
conditions  are  treated  and  results  compared  with  those  obtained  by 
using  PEM  code: 

1.  Rectangular  coordinates  with  uniform  pressure  or 
pressure  gradient  at  the  boundaries,  as  shown  in 
Figure  67. 

2.  Rectangular  coordinates  with  linearly  varying 
pressure  at  the  outer  vertical  boundary  and  uniform 
pressure  or  pressure  gradient  at  other  boundaries  as 
shown  in  Figure  68. 

3.  Cylindrical  polar  coordinates  with  uniform  pressure 
or  pressure  gradient  at  the  boundaries  as  shown  in 
Figure  69. 

Closed  form  solutions  for  all  these  cases  of  boundary 
value  problems  are  obtained  by  Fourier  transform  techniques  and 
are  given  on  page  120  onward.  Numerical  results  for  these 
cases  are  calculated  for  materials  with  following  properties: 

Fluid  Compressibility  (c)  =  1.13  x  10_5  (PSI)-1 
Fluid  Viscosity  (n)  =  0.24  x  10-6  lb/in2-min. 

Solid  Permeability  (K)  =  0.75  x  10-13  in2 
Solid  Porosity  (4>)  =  0.31 

The  problem  solved  for  case  1  is  the  same  as  described  on  page  32. 

The  numerical  results  obtained  by  using  the  closed  form  solution 

and  the  PEM  code  are  shown  in  Figures  70-72.  Figure  70  depicts 

the  variation  of  pressure  along  the  height  of  the  billet  (Figure 

3)  for  x  =  0  and  x  =  1.75".  The  PEM  results  are  plotted  for 

three  values  of  time  step  length,  AT  =  1 . ,  .25  and  0.1,  at  the 

end  of  ten  minutes  impregnation .  This  figure  shows  that  the 

accuracy  of  the  results  by  PEM  code  improves  considerably  with 

the  refinement  of  the  time  step  AT.  In  the  calculation  of  closed  form 


solution  results,  the  summation  of  series  on  page  120  was  done  for 
m=n=25.  The  magnitude  of  the  pressure  P  is  almost  independent  of 
m  and/or  n  beyond  25.  There  exists  an  extremely  good  agreement 
between  the  closed  form  solution  results  and  the  PEM  results  with 
AT  =  0.1  minute.  Figure  71  shows  a  comparison  between  the  PEM 
results  and  the  closed  form  solution  results  at  T  =  25  minutes. 

Here  too,  the  PEM  results  for  AT  =  0.1  minute  are  quite  accurate. 
Figure  72  shows  the  variation  of  pressure,  on  the  basis  of  closed 
form  solution,  at  the  point  (0,0),  of  minimum  pressure  in  the 
billet  with  time.  It  takes  about  200  minutes  before  the  pressure 
at  all  points  of  the  billet  reaches  the  value  of  p  =  89.998. 

Similar  observations  are  made  in  case  2  in  which  a  linearly 
varying  horizontal  pressure  at  the  outer  vertical  boundary  of  the 
billet  is  considered.  Figure  68.  A  comparison  between  PEM  and 
closed  form  results  is  given  in  Figures  73,  74,  and  75.  Table  11 
shows  the  PEM  input  data  for  impregnation  analysis  of  a  billet 
with  a  linearly  varying  horizontal  pressure  at  the  vertical 
boundary  as  shown  in  Figure  68.  The  billet  geometry  considered 
here  is  the  same  as  that  in  the  foregoing  illustration,  case  (1), 

PQ  =  90  and  y  =  45.  Figure  75  shows  the  variation  of  pressure  P  at 
(0,0)  with  time,  calculated  by  the  closed  form  solution.  Thus, 
it  takes  about  200  minutes  before  the  pressure  at  (0,0)  reaches 
the  value  of  P  =  134.997,  close  to  the  applied  pressure  P  =  135 
at  (2,0)  . 

The  problem  solved  for  case  3,  a  billet  in  cylindrical 
polar  coordinates,  is  shown  in  Figure  69.  The  finite  element  grid 
in  the  rz-plane  used  for  the  process  environment  model  is  shown 
in  Figure  76.  The  PEM  input  data  for  this  problem  is  given  in 
Table  12.  The  material  properties  used  for  this  case  are  the  same 
as  those  used  for  the  case  1.  The  closed  form  solution  for  this 
problem  with  relevant  boundary  conditions  is  given  on  page  124.  This 
solution  is  expressed  in  terms  of  Bessel  functions.  For  obtaining 
the  roots  of  equation  4  page  124,  a  CDC  mathematics  library  sub¬ 
routine  was  used.  During  the  computation  of  results,  the  series 
summation  was  done  for  i=m=25.  As  before,  the  results  do  not 
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vary  with  the  increase  in  the  summation  of  number  of  terms  over  m 
or  i  beyond  25.  Figures  77-80  show  the  variation  of  pressure  P 
calculated  by  the  PEM  code  and  the  closed  form  solution  at  different 
points  of  the  billet.  Figure  77  shows  the  effect  of  the  change 
in  the  time  march  step  AT.  Further,  it  shows  that  the  results 
for  time  step  length  AT  =  0.1  minute  are  very  close  to  the  closed 
form  solution  results.  Figures  78  and  79  demonstrate  the  difference 
between  the  PEM  results  with  AT  =  1.0  and  .25  and  the  closed  form 
solution  results  at  T  =  20  and  40  minutes.  Figure  80  gives  the 
change  in  pressure  at  the  point  of  minimum  pressure  with  respect 
to  the  time  T  calculated  on  the  basis  of  the  closed  form  solution. 

It  takes  about  315  minutes  to  reach  a  value  of  P  =  89.636. 


The  governing  equations  for  the  flow  of  compressible  fluid 
through  a  porous  medium  are  given  below: 

The  flux  density  v^.  =  -  (£  =  x,  y,  z) 


where 


»-p0*c(p-po) 


pQ  =  Fluid  density  at  Reference  pressure  PQ 

1  3o 

C  =  Compressibility  constant  =  -  ^ 

K  =  Permeability  of  the  material 
M  =  Fluid  viscosity 

The  flow  of  compressible  liquid  through  an  isotropic  incompressible 
medium  in  terms  of  pressure  is  governed  by  (Darcy's  Law) 


n2P  3  2P  32P  <t>yC  3P 

-  2  -2-2  K  3T 

3x  3y  3z 

(x,y,z)  Cartesian  coordinates,  T 
4>  =  porosity  of  the  material 


time 
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J 


P(x,b,T)  =  P0 


=  0b  P(x,z,o)  =0  P(a,z,-n  =  Pc 


°1 


dP\  =0 

<*z  2=0 


Figure  67.  Coordinate  Axis  and  Pressure  Boundary  Conditions 
for  a  Rectangular  Billet  Impregnation  Analysis. 


Case  1:  Consider  a  two-dimensional  case  in  which  a 
quadrant  of  the  billet  has  been  modeled.  The  initial  and 
boundary  conditions  are  shown  in  Figure  67.  Using  Fourier  trans¬ 
form  technique  the  solution  to  this  problem  is: 


p  _  p  f\  4  \  Y  (~Dn  -a2T/k  V  (-I)™  -(32T/k 

p  -  p0  1_ib  ^  — u —  e  cos  ax  — g —  e  cosBz 

L  n=0  m=0  p 


( 2n+l) 


(2m+l) 


a  *  — n — -  8  - 


z  u  -  H’C 

f  K.  jr 
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Case  2:  Consider  a  2D  case  in  which  the  pressure  boundary 
condition  at  x=a  is  given  by  equation  (1)  as  shown  in  Figure  68. 

P(x,b,T)  =  PA 


P{a,z,T)  = 


PQ  +  Y  (b-z) 


-(1) 


1Z  I 
9x  1 x=0 


92  1 z=0 


P(x,z,0)  =  0 


The  solution  of  this  boundary  value  problem  is: 


2.  2, 


P  ( x ,  z  ,  T)  =  PQ+Y(b-z)-^ 


\ v  \  .  -(A  +y  ) T/k 

>  >  A  CosA  x  Cosy  z  e  m  n 

/  -j  ran  m  n 


m=o  n=o 


where 


A  _  2m+l  _  2n+l 

A  =  — x -  Tr ,  y  =  — jr-: -  7T 

m  2a  n  2b 


A  =  (~Dn  +  _J_  j  ('I)1" 

mn  j  0  y  21  X 

I  n  y  >  m 

n 
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Case  3:  The  field  equation  for  a  cylindrical  billet  in 
Figure  72  can  be  written  as 

13,  3P,  32P  ,  3P  ,, 

r  3?  (r  3F>  +  7T  '  k  3T  ” 

a  z 


The  solution  for  the  boundary  value  problem  given  in  Figure  72  is: 


P(r,z,t)  =  q. 


,  Z  S  (_1)  4b  VbXi,R0(Xir)cos|V  -U?+u2)T/k 

1  pi,  \  o  'o  o  o  e  z  m 

i=l  m=0  u mAi  {b^RQ  ^ ( X.b) -a^Rg (X.^) } 


2m+l 

P  =  — sr> —  TT 
25. 


where 


J  (X  b) 

Rn(Xir)  =  Jn(X.r)  -  Yn(X.r) 


n=0 , 1  (3) 


A  prime  means  the  differentiation  with  respect  to  the  argument 
and  X^  are  real  roots  of: 


J0(Xib)YQ(Xia)  -  JQ  (Xia)YQ(Xib)  =  0. 
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PEM  INPUT  DATA  FOR  IMPREGNATION  ANALYSIS  OF  A  RECTANGULAR  BILLET 
WITH  VARIABLE  HORIZONTAL  PRESSURE 


Figure  75. 


Pressure  at  the  Point  (0,0)  in  the  Billet  Versus 
Time  During  Impregnation,  Cartesan  Coordinates. 
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un.'  78.  Comparison  of  Pressure  Results,  for  Impregnation 
Studies  in  Cylindrical  Polar  Coordinates,  Between 
the  Closed  Form  and  PEM  Calculations  after  Twenty 
M  i  n  u  t  o  s  . 
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Figure  79.  Comparison  of  Results,  for  Impregnation  Studies  in 

Cylindrical  Polar  Coordinates,  Between  the  Closed  Form 
and  PEM  Calculations,  T  =  40  min. 
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Fiqure  80.  Pressure  Variation  at  A  (3.1,  0)  in  the  Billet 
Versus  Time  During  Impregnation  in  Cylindrical 
Coordinates . 
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